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Abstract 

Because  pulsar  phase,  spin  rate,  and  spin-down  rate  cannot  be  determined  a  priori,  a  pulsar’s 
possible  contribution  to  a  time  scale  is  best  measured  using  a  statistic  insensitive  to  these  quantities. 
One  such  statistic  isat(r)  (Matsakis,  Taylor,  and  Eubanks  1996,  hereafter  [1]),  which  is  based  upon 
the  third-order  polynomial  variations  of  the  clock  phases  and  named  in  analogy  to  <t¥(t).  It  is  similar 
to  measures  based  upon  the  variance  of  the  third  differences,  but  more  robust  in  the  presence  of 
data  gaps.  About  6  years  of  data,  from  MJD  47752  to  50343  (18  Aug  89-17  Sep  96),  were  used  to 
determine  the  <t,(t)  of  cesium  and  maser  standards  in  a  variety  of  imperfect  but  probably  acceptable 
ways;  directions  for  improved  reductions  are  also  indicated.  Because  of  the  considerably  improved 
stability  of  cavity -tuned  masers  and  Hewlett-Packard  Model  5071  cesium  standards,  the  best  results 
use  only  the  data  taken  since  MJD  49400  (17  Feb  94),  by  which  time  most  of  the  new  standards 
had  been  acquired.  An  upper  limit  to  the  time  variation  of  the  fine  structure  constant  is  also  given. 


INTRODUCTION 

The  possibility  of  using  the  exceptional  rotational  stability  of  millisecond  pulsars  to  generate 
a  time  scale  has  long  been  of  interest.  I2-1!  Although  the  phase,  rate,  and  spin-down  of  their 
observed  radio  pulsations  are  physically  meaningful  (see  [8]),  they  are  not  predictable  on  the 
basis  of  measurable  quantities,  so  that  their  contribution  to  any  state-of-the-art  terrestrial  time 
scale  can  only  be  to  constrain  derivatives  of  higher  order  than  the  frequency  drift.  While  this 
precludes  their  use  as  absolute  frequency  standards,  incorporating  pulsar  data  in  the  time  scale 
is  still  a  highly  worthy  goal. 

Clock  data  are  commonly  analyzed  using  a  statistic  called  ov,  the  square  root  of  the  “Allan 
variance,”  which  is  proportional  to  the  mean  square  of  the  second  differences  of  a  series  of 
clock  offset  measurements.  Second  differences  are  used  because  they  are  insensitive  to  the 
clock  phase  offset  and  frequency  bias,  and  because  the  frequency  noise  of  these  standards  is 
often  white  over  time  scales  of  less  than  a  few  days  or  weeks.  In  the  presence  of  arbitrary 
frequency  drifts  and  the  “redder”  noise  components  which  dominate  these  time  series  on  longer 
time  scales,  it  is  a  natural  extension  to  consider  statistics  such  as  the  Hadamard  variance,  which 
is  based  upon  the  third  differences.  f9l  Due  to  the  irregular  spacing  of  pulsar  data,  Taylor  I4  J 
suggested  use  of  a  statistic  analogous  to  the  Modified  Hadamard  Variance,  which  he  termed  a „ 
which  is  related  to  a  third-order  fitted  polynomial  of  the  timing  residuals;  the  definition  of  ax 
was  slightly  altered  by  Matsakis,  Taylor,  and  Eubanks.dl  This  statistic  describes  the  lowest-order 
deviations  remaining  in  a  pulsar  time  series  after  the  phase,  frequency,  spin-down  rate,  and 
astrometric  parameters  have  been  determined  by  comparison  with  terrestrial  time,  and  their 
effects  removed. 
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Since  a,  is  inherently  insensitive  to  the  pulsar  properties  which  are  observationally  arbitrary,  it  is 
ideally  suited  for  comparing  pulsar  stabilities  with  those  of  other  time  scales.  Such  comparisons 
are  made  in  [1]  and  elsewhere;  this  paper  presents  an  analysis  of  the  az  characteristics  of  the 
individual  frequency  standards  kept  at  the  U.S.  Naval  Observatory  (USNO). 


THE  USNO  CLOCK  DATA 

The  USNO  clock  data  were  taken  using  the  Digital  Acquisition  System  (DAS),  and  consist 
of  hourly  timing  differences  between  the  “zero  crossings”  of  the  5  MHz  signal  output  by 
individual  frequency  standards  and  the  USNO  Master  Clock  #2  (MC2),  which  is  referred  to 
as  UTC(USNO)  by  the  International  Bureau  of  Weights  and  Measures  (BIPM),  and  in  this 
time  period  was  the  output  of  a  maser  steered  daily  toward  the  weighted  mean  of  the  USNO 
ensemble,  which  itself  is  steered  toward  TAI.  Data  from  70  Hewlett-Packard  Model  5071  cesium 
standards  and  12  cavity-tuned  Sigma-Tau  masers  were  used;  data  from  older-style  frequency 
standards  were  ignored.  The  data  measurement  precision  was  less  than  35  picoseconds  rms 
(Breakiron,  private  communication). 

Individual  clock  differences  with  MC2  were  converted  to  frequencies  and  edited  for  outliers. 
Four  maser  and  five  cesium  clock  time  series  with  jumps  were  subdivided  into  two  series  each 
so  as  to  avoid  imposing  large  rate  or  drift  jumps  in  the  data.  Data  from  one  cesium  clock 
were  subdivided  into  three  series.  For  four  masers  and  ten  cesium  clocks,  the  need  to  break  a 
time  series  up  was  not  dear,  and  so  both  the  original  series  and  the  “optional”  subsets  were 
retained.  The  number  of  standards  plotted  is,  therefore,  slightly  larger  than  the  number  of 
physical  clocks  we  had,  and  about  double  the  number  of  standards  currently  used  by  the  USNO 
to  compute  its  time  scales. 


PAIRWISE  COMPUTATION  OF  az 

Once  the  editing  was  complete,  the  a,  statistic  was  generated  for  the  difference  series  generated 
by  subtracting  each  clock’s  frequency  data  from  other  clock  data,  from  MC2,  from  the  USNO 
mean,  from  the  free-running  USNO  time  scale  A.1,  from  the  spline-interpolated  MC2-TAI,  and 
from  the  spline-interpolated  TT(96)  time  series,  which  is  generated  by  the  BIPM  particularly 
for  pulsar  work  (see  [10]). 

To  compute  (and  define)  <xz  using  frequency  data  instead  of  timing  residual  data,  the  recipe 
given  in  [10]  was  adapted  as  follows: 

1.  For  each  clock,  phase  data  pairs  separated  by  1  day  were  differenced  to  generate  the 
frequencies  y(tf). 

2.  The  data  were  divided,  according  to  time  of  measurement,  into  subsequences  defined  by 
continuous  intervals  of  length  r. 

3.  By  minimizing  the  weighted  sum  of  squared  differences,  [(y(tf)  -  yfo))]8,  the  data  in  each 
subsequence  were  modelled  to  the  derivative  of  the  cubic  function  used  by  [1]: 

Y(t)  —  c i  +  2  x  pj(t,  -  to)  +  3  x  c3(tj  -  to)s  (1) 

The  weights  were  chosen  to  be  unity.  All  overlapping  subsequences  were  included  if  they 
contained  four  or  more  measurements  and  the  interval  between  first  and  last  measurement  was 


at  least  r/v/2. 

3.  It  follows  that 


(2) 


1/2 


where  angle  brackets  denote  averaging  over  the  subsequences,  weighted  by  the  inverse  squares 
of  the  formal  errors  in  c3. 

A  comparison  with  analyses  based  upon  the  Hadamard  variance,  such  as  Hutsellt11!,  Hutsell  et 
aU12l,  and  Walter!13-14!,  can  be  made  by  noting  that  the  “square  root  of  the  modified  Hadamard 
variance”  has  a  r-dependence  similar  to  that  of  nx,  after  allowance  is  made  for  the  factor  of 
s/3/10  normalization  difference  and  for  the  fact  that  the  r  used  here  is  three  times  larger  than 
the  t  used  to  define  the  Hadamard  variance. 


ESTIMATION  OF  INDIVIDUAL  AND  AVERAGE  az 

If  one  is  willing  to  make  the  questionable  assumption,  which  is  discussed  in  the  next  section, 
that  the  clock  errors  are  independent  and  uncorrelated,  then  the  associated  with  each  clock 
difference  is  equal  to  the  sum  of  the  “intrinsic”  a\  variances  of  each  clock.  It  then  becomes 
mathematically  possible  to  derive  the  individual  variances  from  a  least-squares  fit  to  all  pairs 
via  an  N-comered-hat  analysis,  which  is  a  weighted  least-squares  solution  to  find  those  values 
for  individual  clock  crj’s  which  best  fit  the  statistics  of  the  clock  pair-difference  <r®’s. 

In  the  solutions  individual  clock  pairs  can  be  weighted  in  several  ways.  Solutions  were  generated 
by  weighting  each  pair  identically,  by  using  the  statistical  weights  derived  from  variances  near 
r  of  a  month  when  compared  to  MC2,  by  downweighting  or  deweighting  masers  or  cesiums, 
by  excluding  all  but  the  best  third  of  any  type  (as  measured  by  their  deviations  from  MC2  at 
intermediate  periods),  and  by  combinations  of  these.  While  there  was  not  a  large  difference 
between  the  results  of  these  weighting  schemes  for  any  given  clock,  the  lowest  individual 
variances  usually  were  those  inferred  from  the  formal  errors  generated  by  assuming  an  identical 
“white”  error  dominating  all  measurements,  so  that  the  results  depend  only  on  the  irregularities 
in  data  spacing. 

Different  solutions  were  found  using  subsets  of  the  data.  Lower  variances  were  found  if  the 
data  were  restricted  to  be  after  MJD  49400,  and  still  lower  were  found  by  excluding  clocks 
which  did  not  completely  fill  the  interval  from  MJD  49400  to  50343. 

In  Figures  1  a-d  the  data  for  the  individual  cesium  standards  are  plotted,  using  only  data  after 
MJD  49400  and  weighting  each  clock  pair  by  the  formal  error  with  MC2.  Figure  2  is  a  similar 
plot  for  the  masers.  Note  that  individual  standards  can  differ  by  an  order  of  magnitude. 

In  Figures  3  and  4,  various  averages  of  the  (non-imaginary)  derived  values  of  oz  are  plotted. 
Most  of  the  curves  are  from  estimates  derived  with  different  weighting  schemes.  Denoted 
with  a  “1”  is  the  average  when  the  standards  are  compared  to  MC2,  while  that  with  a  “2” 
is  the  average  comparing  the  standards  to  TAI.  The  curves  are  presented  as  an  indicator  of 
robustness.  The  limits  to  the  technique  are  shown  by  the  fact  that  the  curves  in  Figure  3, 
which  incorporate  all  the  data,  are  higher  than  those  in  Figure  4,  which  are  based  on  only  data 
taken  more  recently  than  MJD  49400.  The  estimates  in  Figures  3  and  4  are  not  corrected  for 
the  bias  due  to  the  fact  that  the  average  measured  value  of  <r2  has  a  chi-squared  distribution 
and  for  this  distribution  the  median-centered  logarithm  of  <rz  is  lower  than  the  average  value 
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by  .017/n,  where  n  is  the  number  of  degrees  of  freedom. ril  Applying  this  correction  yields 
median-centered  curves  lower  than  the  plotted  figures.  The  difference,  in  logarithmic  units,  is 
.17  for  the  longest  r  (7.9),  and  half  as  much  for  r  =  7.6. 

ROBUSTNESS  OF  THE  SOLUTIONS 

The  N-cornered-hat  method  will  not  reveal  any  “external”  error  which  is  correlated  with  the 
ensemble  average,  and  it  is  sensitive  to  the  fact  that  the  observed  variances  are  themselves 
random  variables  (see  [15]).  It  also  can  give  highly  erroneous  results  in  the  presence  of  large 
measurement  noise,  data  gaps,  nonstationarity,  and  nonzero  correlations  between  clock  outputs. 
Such  problems  often  lead  to  the  derivation  of  negative  estimates  for  individual  clock  <r3s  (which 
were  discarded  in  this  work).  The  derivation  of  negative  variances  increases  with  r  because 
the  clock  correlations  fail  to  average  to  zero  as  the  number  of  independent  points  decreases. 
Some  benefit  is  gained  from  the  stochastic  character  added  to  the  covariance  terms  by  data 
gaps  and  other  irregularities. 

In  any  event,  the  solution  results  do  not  appear  to  vary  significantly  when  different  weighting 
functions  are  used  in  the  fit,  and  they  are  consistent  with  comparisons  between  each  clock  and 
TAI  or  the  various  USNO  time  series  which  are  more  accurate  and  stable  than  their  components. 
Not  surprisingly,  the  scatter  between  the  different  methods  increases  at  the  longest  intervals, 
which  is  the  regime  where  the  approximations  are  most  questionable. 

It  is  possible  to  generate  improved  clock  variances  by  making  a  definable  and  explorable  set 
of  assumptions  about  the  characteristics  of  just  one  reference  clock!16-19!  or  through  use  of  a 
reference  time  scale  much  more  precise  than  the  individual  clocks.  Using  these  variances,  it  is 
possible  to  estimate  the  characteristics  of  the  noise.  Although  the  presence  of  non-white  noise 
is  obvious  in  all  the  figures,  we  have  not  attempted  to  quantify  its  characteristics.  This  is  best 
done  by  a  multivariance  approach  I13l,  which  is  based  upon  the  use  of  complementary  statistical 
measures,  such  as  <rv.  We  intend  to  address  these  problems  in  subsequent  papers. 


AN  ASIDE  ON  COSMOLOGICAL  IMPLICATIONS 

It  has  been  pointed  outP°l  that  variations  in  the  fine  structure  constant,  such  as  those  predicted 
by  recent  grand  unified  theories  and  Dirac’s  Large  Numbers  Hypothesis,  would  lead  to  a 
differential  drift  in  frequency  standards  which  are  based  upon  different  atomic  transitions.  The 
frequency  drifts  of  best-measured  individual  standards  in  this  analysis  range  over  3  xlO~13/year 
peak  to  peak,  and  we  find  no  significant  difference  between  the  two  types  of  standards  to  a 
1-sigma  upper  limit  of  2.5xl0_14/year.  According  to  the  formulas  of  [20],  this  corresponds  to  a 
limit  of  3.2xl0~14/year  in  the  fractional  variation  of  the  fine  structure  constant,  and  is  similar  to  a 
limit  they  derive  by  comparing  mercury-  and  hydrogen-based  standards.  It  is  possible  that  a  large 
variation  in  the  fine  structure  constant  is  being  fortuitously  masked  in  our  data  by  some  equally 
large  (differential)  instrumental  effect  common  to  one  type  of  standard  as  maintained  at  the 
USNO.  We  consider  this  to  be  extremely  unlikely  because,  if  it  were  true,  a  different  differential 
instrumental  effect  (of  a  specific  magnitude)  would  also  have  to  be  present  in  the  data  analyzed 
by  [20],  which  involved  different  atomic  elements  and  different  institutions-the  mercury-based 
clocks  were  built  by  and  maintained  at  the  Jet  Propulsion  Laboratory  and  the  primary  cesium 
standards  were  built  and  maintained  by  the  Physikalisch-Technische  Bundesanstalt  (PTB)J211 
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COMPARISONS  WITH  PULSAR  DATA 


The  purpose  of  this  paper  is  to  present  the  observed  clock  statistics,  but  no  discussion  of  those 
statistics  would  be  complete  without  a  comparison  with  pulsar  data.  In  the  last  figure,  data 
from  Figures  3  and  4  of  [l]  are  combined.  The  az  statistic  of  the  timing  difference  between 
two  pulsars  and  TAI  are  presented,  as  is  the  comparison  between  the  free-running  terrestrial 
time  scales  of  the  USNO  (A.l)  and  the  PTB,  to  which  the  BIPM  steered  TAI  over  the  time 
range  of  the  pulsar  data.  To  generate  these  curves,  data  from  MJD  44979  to  50079  were  used 
(more  recent  data  were  not  used  because  no  high-quality  pulsar  data  are  available  to  us  due 
to  the  Arecibo  telescope  ’  upgrade  work).  Although  one  would  infer  from  those  curves  that 
pulsar  data  can  contribute  meaningfully  to  time  scales  over  periods  exceeding  a  few  years,  the 
recent  incorporation  of  HP  5071  cesiums  and  cavity-tuned  masers  into  the  terrestrial  time  scales 
results  in  none  of  the  comparisons  being  stationary  in  a  statistical  sense.  The  lowest-lying  curve 
shows  the  difference  between  the  two  terrestrial  time  scales  when  only  the  last  2  years  of  data 
are  used,  which  corresponds  to  the  introduction  of  the  new  standards.  Since  the  pulsar  data 
analysis  necessarily  depends  upon  comparison  with  the  older  time  scales,  an  adequate  pulsar 
database  for  comparison  will  not  be  available  for  several  more  years. 
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Figure  1  a.  Twenty  individual  cesium  5071s,  MJD  49400-50343 
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Figure  1  b.  Twenty  more  individual  cesium  5071s,  MJD  49400-50343 
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Figure  3  b.  Averages  of  Sigma-Tau  masers,  MJD  47752-50343 


Figure  4  b.  Averages  of  Sigma-Tau  masers,  MJD  49400-50343 


Figure  5.  Pulsar  and  terrestrial  time  scale  stabilities 
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